5  Chapter 3 Exercise Solutions

5.1 Solution to Exercise 3.1

Label 365 cards with the numbers 1 through 365. Shuffle the cards and deal \(n\) with replacement. If the \(n\) cards dealt all have different numbers, then \(B\) does not occur; otherwise \(B\) occurs.

5.2 Solution to Exercise 3.2

Create a spinner with two sectors, one sector labeled “make” which accounts for 40% of the area, and the other 60% labeled “miss”. Spinner the spinner until it lands on “make” and then stop; let \(X\) be the total number of spins.

5.3 Solution to Exercise 3.3

Create a spinner with three equal sectors, labeled 1, 2, 3. Spin the spinner 3 times. Let \(X\) be the number of distinct numbers the spinner lands on; for example if the spinner lands on 3 then 1 then 3, then \(X\) is 2. Let \(Y\) be the number of spins that land on 1.

5.4 Solution to Exercise 3.4

  1. A Uniform(0, 12) spinner is basically a clock. It goes from 0 to 12 in equally spaced increments. See Figure 5.1.
  2. To simulate \(X\), flip the coin to determine if the point is in the “east” or “west” and then spin the Uniform(0, 12) spinner to determine how far from the center in the east or west direction. Let the result of the spin be \(U_1\). If the coin lands on heads (east) set \(X=U_1\); if the coin lands on tails (west) set \(X=-U_1\). To simulate \(Y\), spin the spinner again to get \(U_2\), and flip the coin again. If the result of the flip is heads (north) then set \(Y = U_2\); if the result of the flip is tails (south) then set \(Y=-U_2\). Compute the distance from \((X, Y)\) to (0, 0) as \(R=\sqrt{X^2+Y^2}\). If \(R>12\) the dart would land off the board, so just discard the repetition and try again.
  3. Exercise 2.13 shows that \(R\) is not uniformly distributed between 0 and 12. In particular, \(R\) is more likely to be between 11 and 12 than it is to be between 0 and 1. Therefore a Uniform(0, 12) spinner would not be appropriate.
Figure 5.1: A continuous Uniform(0, 12) spinner. Only selected rounded values are displayed, but in the idealized model the spinner is infinitely precise so that any real number between 0 and 12 is a possible outcome.

5.5 Solution to Exercise 3.5

Figure 5.2: Spinner representing the distribution of \(R\) in Exercise 3.5

5.6 Solution to Exercise 3.6

P = BoxModel([1, 2, 3], size = 3)

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))

(RV(P) & X & Y).sim(10)
Index Result
0((1, 2, 1), 2, 2)
1((3, 1, 3), 2, 1)
2((1, 1, 3), 2, 2)
3((2, 2, 3), 2, 0)
4((3, 1, 2), 3, 1)
5((3, 3, 2), 2, 0)
6((2, 1, 1), 2, 2)
7((2, 3, 2), 2, 0)
8((1, 2, 1), 2, 2)
......
9((3, 3, 2), 2, 0)

5.7 Solution to Exercise 3.7

P = BoxModel([1, 2, 3], size = 3, probs = [0.1, 0.3, 0.6])

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))

(RV(P) & X & Y).sim(10)
Index Result
0((1, 2, 3), 3, 1)
1((2, 3, 3), 2, 0)
2((3, 1, 2), 3, 1)
3((2, 3, 3), 2, 0)
4((3, 3, 3), 1, 0)
5((1, 2, 3), 3, 1)
6((3, 3, 2), 2, 0)
7((3, 2, 3), 2, 0)
8((3, 2, 1), 3, 1)
......
9((2, 3, 2), 2, 0)

5.8 Solution to Exercise 3.8

The code below defines a random variable \(X\) that counts the number of distinct birthdays in the group of \(n\). Event \(B\) occurs if \(X<n\).

n = 30

P = BoxModel(list(range(365)), size = n)

def count_distinct(u):
  return len(set(u))

X = RV(P, count_distinct)

B = (X < n)

B.sim(10)
Index Result
0True
1True
2True
3False
4False
5True
6False
7True
8False
......
9False

5.9 Solution to Exercise 3.9

p_success = 0.4

# define a function that takes as an input a sequence of 0s and 1s (omega)
# and returns when the first 1 occurs
def count_until_first_success(omega):
    for i, w in enumerate(omega):
        if w == 1:
            return i + 1 # the +1 is for zero-based indexing

# Either of the following probability spaces works        
# P = Bernoulli(p_success) ** inf
P = BoxModel([1, 0], probs = [p_success, 1 - p_success], size = inf)

X = RV(P, count_until_first_success)

(RV(P) & X).sim(10)
Index Result
0((0, 1, 0, 0, 1, 0, ...), 2)
1((0, 1, 1, 0, 0, 0, ...), 2)
2((0, 0, 0, 0, 1, 1, ...), 5)
3((0, 0, 0, 1, 0, 0, ...), 4)
4((0, 0, 0, 1, 0, 1, ...), 4)
5((0, 0, 0, 0, 1, 0, ...), 5)
6((0, 0, 1, 1, 0, 1, ...), 3)
7((1, 0, 0, 0, 0, 1, ...), 1)
8((1, 1, 1, 1, 0, 0, ...), 1)
......
9((0, 1, 0, 1, 0, 1, ...), 2)

5.10 Solution to Exercise 3.10

In the code below

  • Uniform(0, 12) ** 2 corresponds to the two spins
  • BoxModel([-1, 1], size = 2) corresponds to the two coin flips
  • Uniform(0, 12) ** 2 * BoxModel([-1, 1], size = 2) simulates a pair which we unpack and define as RVs U, F, but U is itself a pair and so is F.
  • U[0] is the first spin and F[0] is the first flip.
U, F = RV(Uniform(0, 12) ** 2 * BoxModel([-1, 1], size = 2))

X = U[0] * F[0]

Y = U[1] * F[1]

R = sqrt(X ** 2 + Y ** 2)

(U & F & X & Y & R).sim(10)
Index Result
0((0.9251290000106569, 3.643053509029234), (-1, 1), -0.9251290000106569, 3.643053509029234, 3.7586836...
1((0.9933298365253349, 11.677550742000548), (-1, -1), -0.9933298365253349, -11.677550742000548, 11.71...
2((8.60453131409852, 7.019557104396412), (1, 1), 8.60453131409852, 7.019557104396412, 11.104599996271...
3((1.984668527085105, 5.259208171690422), (-1, -1), -1.984668527085105, -5.259208171690422, 5.6212258...
4((7.75428773226977, 2.056109940529828), (1, 1), 7.75428773226977, 2.056109940529828, 8.0222544413883...
5((3.9392050687422318, 9.170546367635616), (1, -1), 3.9392050687422318, -9.170546367635616, 9.9807944...
6((4.213564925337111, 1.949138900670969), (-1, -1), -4.213564925337111, -1.949138900670969, 4.6425501...
7((11.844683022909935, 7.334515271868463), (-1, -1), -11.844683022909935, -7.334515271868463, 13.9316...
8((2.6764651254764407, 3.851646367314334), (1, -1), 2.6764651254764407, -3.851646367314334, 4.6902713...
......
9((6.791001181714114, 11.04037826877519), (1, 1), 6.791001181714114, 11.04037826877519, 12.9617764742...

Currently the \((X, Y)\) pairs are uniformly distributed in the box with sides [-12, 12]. We’ll see how to discard points off the board later.

plt.figure()
(X & Y).sim(1000).plot()
plt.show()

5.11 Solution to Exercise 3.11

P = Uniform(1, 4) ** 2

X = RV(P, sum)
Y = RV(P, max)

(RV(P) & X & Y).sim(10)
Index Result
0((1.2782363749769647, 1.064015488396242), 2.342251863373207, 1.2782363749769647)
1((2.4126689042489167, 2.779910314685748), 5.1925792189346645, 2.779910314685748)
2((1.709716471619461, 1.5353605853331378), 3.245077056952599, 1.709716471619461)
3((2.334942785988934, 1.238657119642979), 3.5735999056319128, 2.334942785988934)
4((1.4440431323205587, 1.5631192654444537), 3.007162397765012, 1.5631192654444537)
5((3.8552598615279177, 1.3649673196981444), 5.220227181226062, 3.8552598615279177)
6((3.859363812790382, 3.300649080165231), 7.1600128929556135, 3.859363812790382)
7((3.897663679605871, 1.7288994707957444), 5.626563150401616, 3.897663679605871)
8((3.1238652074720403, 3.3460330553428697), 6.46989826281491, 3.3460330553428697)
......
9((1.2332935082600125, 2.685253288839661), 3.9185467970996735, 2.685253288839661)
plt.figure()
(X & Y).sim(1000).plot()
plt.show()

5.12 Solution to Exercise 3.12

Part 1

See previous solutions which discuss how to simulate \((X, Y)\) pairs.

Simulate many \((X, Y)\) pairs, say 10000.

  1. To approximate \(\textrm{P}(X = 2)\): Divide the number of repetitions with \(X = 2\) by the total number of repetitions
  2. To approximate \(\textrm{P}(Y = 1)\): Divide the number of repetitions with \(Y = 1\) by the total number of repetitions
  3. To approximate \(\textrm{P}(X = 2, Y = 1)\): Divide the number of repetitions with \(X = 2\) and \(Y=1\) by the total number of repetitions

Part 2

P = BoxModel([1, 2, 3], size = 3)

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
x.count_eq(2) / x.count()
0.6801
y.count_eq(1) / x.count()
0.4453
((x == 2) * (y == 1)).sum() / x.count()
0.2275

Part 3

x.tabulate(normalize = True)
Value Relative Frequency
10.1021
20.6801
30.2178
Total1.0
plt.figure()
x.plot()
plt.show()

y.tabulate(normalize = True)
Value Relative Frequency
00.2877
10.4453
20.2314
30.0356
Total1.0
plt.figure()
y.plot()
plt.show()

x_and_y.tabulate(normalize = True)
Value Relative Frequency
(1, 0)0.0665
(1, 3)0.0356
(2, 0)0.2212
(2, 1)0.2275
(2, 2)0.2314
(3, 1)0.2178
Total1.0
plt.figure()
x_and_y.plot('tile')
plt.show()

5.13 Solution to Exercise 3.13

Part 1

See previous solutions which discuss how to simulate \((X, Y)\) pairs.

Simulate many \((X, Y)\) pairs, say 10000.

  1. To approximate \(\textrm{P}(X = 2)\): Divide the number of repetitions with \(X = 2\) by the total number of repetitions
  2. To approximate \(\textrm{P}(Y = 1)\): Divide the number of repetitions with \(Y = 1\) by the total number of repetitions
  3. To approximate \(\textrm{P}(X = 2, Y = 1)\): Divide the number of repetitions with \(X = 2\) and \(Y=1\) by the total number of repetitions

Part 2

P = BoxModel([1, 2, 3], size = 3, probs = [0.1, 0.3, 0.6])

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
x.count_eq(2) / x.count()
0.6534
y.count_eq(1) / x.count()
0.2371
((x == 2) * (y == 1)).sum() / x.count()
0.1323

Part 3

x.tabulate(normalize = True)
Value Relative Frequency
10.2418
20.6534
30.1048
Total1.0
plt.figure()
x.plot()
plt.show()

y.tabulate(normalize = True)
Value Relative Frequency
00.7356
10.2371
20.026
30.0013
Total1.0
plt.figure()
y.plot()
plt.show()

x_and_y.tabulate(normalize = True)
Value Relative Frequency
(1, 0)0.2405
(1, 3)0.0013
(2, 0)0.4951
(2, 1)0.1323
(2, 2)0.026
(3, 1)0.1048
Total1.0
plt.figure()
x_and_y.plot('tile')
plt.show()

5.14 Solution to Exercise 3.14

Part 1

See previous solutions for description of how to simulate values of \(X\).

Simulate many values of \(X\) and summarize the simulate values and their relative frequencies to approximate the distribution of \(X\).

To approximate \(\textrm{P}(X > 3)\): divide the number of repetitions with \(X>3\) by the total number of repetitions.

Part 2

p_success = 0.4

# define a function that takes as an input a sequence of 0s and 1s (omega)
# and returns when the first 1 occurs
def count_until_first_success(omega):
    for i, w in enumerate(omega):
        if w == 1:
            return i + 1 # the +1 is for zero-based indexing

# Either of the following probability spaces works        
# P = Bernoulli(p_success) ** inf
P = BoxModel([1, 0], probs = [p_success, 1 - p_success], size = inf)

X = RV(P, count_until_first_success)

x = X.sim(10000)

x
Index Result
01
12
21
32
45
52
61
71
81
......
99995
x.tabulate(normalize = True)
Value Relative Frequency
10.3962
20.2409
30.1433
40.0858
50.053
60.0347
70.0184
80.0124
90.0072
100.003
110.0019
120.0016
130.0005
140.0006
150.0002
160.0001
170.0001
210.0001
Total1.0
plt.figure()
x.plot()
plt.show()

x.count_gt(3) / x.count()
0.2196

5.15 Solution to Exercise 3.15

Part 1

See previous solutions for how to simulate \((X, Y)\) pairs.

Simulate many \((X, Y)\) pairs, say 10000. To approximate:

  1. \(\textrm{P}(X < 3.5)\): divide number of repetitions with \(X<3.5\) by the total number of repetitions
  2. \(\textrm{P}(Y > 2.7)\): divide number of repetitions with \(Y>2.7\) by the total number of repetitions
  3. \(\textrm{P}(X < 3.5, Y > 2.7)\): divide number of repetitions with \(X<3.5\) and \(Y>2.7\) by the total number of repetitions

Part 2

We usually summarize simulated values of continuous random variables with histogram.

P = Uniform(1, 4) ** 2

X = RV(P, sum)
Y = RV(P, max)

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
plt.figure()
x.plot()
plt.show()

plt.figure()
y.plot()
plt.show()

plt.figure()
x_and_y.plot()
plt.show()

plt.figure()
x_and_y.plot('hist')
plt.show()

5.16 Solution to Exercise 3.16

Simulate 10000 repetitions and find the relative frequency of repetitions on which at least 2 people share a birthday. The margin of error is \(1/\sqrt{10000} = 0.01\)

n = 30

P = BoxModel(list(range(365)), size = n)

def count_distinct(u):
  return len(set(u))

X = RV(P, count_distinct)

B = (X < n)

p_hat = B.sim(10000).count_eq(True) / 10000

p_hat
0.704

We estimate (with 95% confidence) that the probability that at least two people in a group of 30 share the same birthday is 0.704 with a margin of error of 0.01. In other words, we estimate (with 95% confidence) that the probability that at least two people in a group of 30 share the same birthday is between 0.69 and 0.71.

5.17 Solution to Exercise 3.17

The margin of error for any single probability would be \(1/\sqrt{10000} = 0.01\), but since we are estimating many probabilities simultaneously we would use a margin of error of \(2\times 0.01\). For example, if the simulated relative frequency of \(X = 2\) is 0.662, then we approximate \(\textrm{P}(X = 2)\) to be 0.662 with a margin of error of 0.02; in other words, we estimate that \(\textrm{P}(X = 2)\) is between 0.642 and 0.682.

5.18 Solution to Exercise 3.18

  1. \(X\) takes the value 210 with probability 1/5 and 35 with probability 4/5.
  2. \(\text{E}(X) = 210(1/5) + 35(4/5)=70\) is the average number of students per class from the instructors’ perspective. If we randomly select an instructor, record the number of students in the instructor’s class, and repeat, the long run average number of students will be 70.
  3. \(\text{P}(X = \text{E}(X)) = \text{P}(X = 70) = 0\). No instructor has a class whose size is equal to the average class size (from the instructor’s perspective).
  4. \(Y\) takes the value 210 with probability 210/350 and 35 with probability 140/350.
  5. \(\text{E}(Y) = 210(210/350) + 35(140/350)=140\) is the average number of students per class from the students’ perspective. If we randomly select a student, record the number of students in the selected student’s class, and repeat, the long run average number of students will be 140.
  6. \(\text{P}(Y = \text{E}(Y)) = \text{P}(Y = 140) = 0\). No student is in a class whose size is equal to the average class size (from the students’ perspective).
  7. Colleges usually report average class size from the instructor/class perspective, which in this case would be 70 students. From the students’ perspective, the average class size is not even close to 70! In fact, it’s twice that size. Some students (140 of them, which is 40% of the total of 350 students) have the benefit of a small class size of 35. But most students (210 of them, which is 60% of the students) are stuck in a large class of 210 students. In other words, most students would be pretty seriously misled if they chose this college based on the advertised average class size of 35 students per class. From the students’ perspective, it seems that 140 is the more relevant average to report. However, neither average really adequately represents what is happening here: there is wide variability in class sizes between large and small.

5.19 Solution to Exercise 3.19

See previous for solutions for how to simulate \((X, Y)\) pairs

Part 1

Simulate many \((X, Y)\) pairs, say 10000

  1. To approximate \(\textrm{E}(X)\): average the simulated \(X\) values (sum the 10000 simulated \(X\) values and divide by 10000)
  2. To approximate \(\textrm{E}(Y)\): average the simulated \(Y\) values (sum the 10000 simulated \(Y\) values and divide by 10000)
  3. To approximate \(\textrm{E}(X^2)\): average the simulated \(X^2\) values (square each simulated \(X\) value then sum the 10000 \(X^2\) values and divide by 10000)
  4. To approximate \(\textrm{E}(XY)\): average the simulated \(XY\) values (for each \((X, Y)\) pair compute the product \(XY\) then sum the 10000 \(XY\) values and divide by 10000)

Part 2

P = BoxModel([1, 2, 3], size = 3)

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
x.mean()
2.118
y.mean()
0.9945
(x ** 2).mean()
4.8078
(x * y).mean()
2.1016

5.20 Solution to Exercise 3.20

See previous for solutions for how to simulate \((X, Y)\) pairs

Part 1

Simulate many \((X, Y)\) pairs, say 10000

  1. To approximate \(\textrm{E}(X)\): average the simulated \(X\) values (sum the 10000 simulated \(X\) values and divide by 10000)
  2. To approximate \(\textrm{E}(Y)\): average the simulated \(Y\) values (sum the 10000 simulated \(Y\) values and divide by 10000)
  3. To approximate \(\textrm{E}(X^2)\): average the simulated \(X^2\) values (square each simulated \(X\) value then sum the 10000 \(X^2\) values and divide by 10000)
  4. To approximate \(\textrm{E}(XY)\): average the simulated \(XY\) values (for each \((X, Y)\) pair compute the product \(XY\) then sum the 10000 \(XY\) values and divide by 10000)

Part 2

P = Uniform(1, 4) ** 2

X = RV(P, sum)

Y = RV(P, max)

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
x.mean()
4.974025031400507
y.mean()
2.9893384679335866
(x ** 2).mean()
26.230541386904225
(x * y).mean()
15.616798267571777

5.21 Solution to Exercise 3.21

See previous parts for how to simulate \(X\) values.

Part 1

To approximate \(\textrm{Var}(X)\):

  1. Simulate many \(X\) values and compute the average of the simulated values
  2. Square each simulated \(X\) value and compute the average of the squared values
  3. Subrtract the square of the average from part a) from the average of the squares in part b).

Part 2

P = BoxModel([1, 2, 3], size = 3)

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

x = X.sim(10000)
(x ** 2).mean() - (x.mean()) ** 2
0.32664375000000057
x.var()
0.32664375000000007
x.sd()
0.5715275583906695

5.22 Solution to Exercise 3.22

Part 1

P = Uniform(1, 4) ** 2

X = RV(P, sum)

Y = RV(P, max)

x_and_y = (X & Y).sim(10000)

x = x_and_y[0]

y = x_and_y[1]
x.plot()

y.plot()

Part 2

\(X\) will have larger standard deviation than \(Y\). The values of \(X\) range from 2 to 8 while the values of \(Y\) only range from 1 to 4. Also, the values of \(Y\) tend to be closer to their mean (around 3.0) than the values of \(X\) are to their mean (around 5).

The values of \(X\) can be anywhere from 0 to 3 units away from the mean of 5, but more values are closer than farther; we might estimate the SD to be less than 1.5. (Remember, the squaring then averaging then taking the square root makes it hard to guess the actual number)

We would estimate the SD of \(Y\) to be less that of \(X\), and maybe a little less than 1, since a large percentage of values are between 2 and 4 (hence less than 1 unit away from the mean).

Part 3

To approximate \(\textrm{Var}(X)\):

  1. Simulate many \(X\) values and compute the average of the simulated values
  2. Square each simulated \(X\) value and compute the average of the squared values
  3. Subrtract the square of the average from part a) from the average of the squares in part b).

Take the square of the variance to approximate the standard deviation.

Similarly for \(Y\).

Part 4

x.var(), x.sd()
(1.4790855015065627, 1.2161765914153104)
y.var(), y.sd()
(0.49498261784722364, 0.7035500109069885)

5.23 Solution to Exercise 3.23

TBA

5.24 Solution to Exercise 3.24

TBA

5.25 Solution to Exercise 3.25

Approximate the conditional distribution of \(Y\) given \(X=1\):

  • Simulate \((X, Y)\) pair
  • If \(X\neq1\) discard the repetition, otherwise keep
  • Repeat until 10000 repetitions with \(X=1\) are obtained
  • Summarize the simulated values of \(Y\) from these repetitions to approximate the conditional distribution of \(Y\) given \(X = 1\)
  • To approximate \(\textrm{P}(Y = 0 | X=1)\) count the number of repetitions with \(Y=0\) and divide by 10000 (remember, all of the repetitions have \(X = 1\))

Approximate the conditional distribution of \(X\) given \(Y=0\) similarly, with a different simulation of 10000 \((X, Y)\) pairs that satisfy \(Y=0\).

P = BoxModel([1, 2, 3], size = 3)

def count_distinct(u):
    return len(set(u))

X = RV(P, count_distinct)

Y = RV(P, count_eq(1))
y_given_Xeq1 = ( Y | (X == 1) ).sim(10000)
y_given_Xeq1.tabulate()
Value Frequency
06727
33273
Total10000
y_given_Xeq1.tabulate(normalize=True)
Value Relative Frequency
00.6727
30.3273
Total1.0
plt.figure()
y_given_Xeq1.plot()
plt.show()

x_given_Yeq0 = ( X | (Y == 0) ).sim(10000)
x_given_Yeq0.tabulate()
Value Frequency
12489
27511
Total10000
x_given_Yeq0.tabulate(normalize=True)
Value Relative Frequency
10.2489
20.7511
Total1.0
plt.figure()
x_given_Yeq0.plot()
plt.show()

5.26 Solution to Exercise 3.26

In the code below

  • Uniform(0, 12) ** 2 corresponds to the two spins
  • BoxModel([-1, 1], size = 2) corresponds to the two coin flips
  • Uniform(0, 12) ** 2 * BoxModel([-1, 1], size = 2) simulates a pair which we unpack and define as RVs U, F, but U is itself a pair and so is F.
  • U[0] is the first spin and F[0] is the first flip.
U, F = RV(Uniform(0, 12) ** 2 * BoxModel([-1, 1], size = 2))

X = U[0] * F[0]

Y = U[1] * F[1]

R = sqrt(X ** 2 + Y ** 2)

( (U & F & X & Y & R) | (R < 12) ).sim(10)
Index Result
0((7.344527801629187, 7.888580359400663), (1, -1), 7.344527801629187, -7.888580359400663, 10.77830175...
1((7.165111597122447, 7.437939755974995), (1, 1), 7.165111597122447, 7.437939755974995, 10.3277186257...
2((7.82900894617325, 2.5244912968145967), (-1, -1), -7.82900894617325, -2.5244912968145967, 8.2259611...
3((2.7441378519691613, 7.577533059576811), (-1, -1), -2.7441378519691613, -7.577533059576811, 8.05911...
4((2.204355482367867, 7.328625099822039), (-1, -1), -2.204355482367867, -7.328625099822039, 7.6529686...
5((8.232590014603588, 7.128198242363699), (1, -1), 8.232590014603588, -7.128198242363699, 10.88975429...
6((6.804211844559606, 6.890857905794813), (-1, 1), -6.804211844559606, 6.890857905794813, 9.684070502...
7((10.488413747074533, 0.9711958594564449), (1, 1), 10.488413747074533, 0.9711958594564449, 10.533282...
8((7.424584555130927, 2.410266921312378), (-1, -1), -7.424584555130927, -2.410266921312378, 7.8060132...
......
9((10.152942855566417, 0.8722243909538392), (-1, -1), -10.152942855566417, -0.8722243909538392, 10.19...
plt.figure()
( (X & Y) | (R < 12) ).sim(1000).plot()
plt.show()

r = ( R | (R < 12) ).sim(10000)

plt.figure()
r.plot()
plt.show()

r.count()
10000
r.count_lt(1) / r.count()
0.0066
r.count_gt(11) / r.count()
0.1649
r.mean()
8.04817722728152
r.sd()
2.825768767689035

5.27 Solution to Exercise 3.27

We condition on \(X\) being “close to” 3.5, say within 0.1 of 3.5

P = Uniform(1, 4) ** 2

X = RV(P, sum)

Y = RV(P, max)
( (X & Y) | (abs(X - 3.5) < 0.1) ).sim(10)
Index Result
0(3.4906064787036755, 2.198778212360655)
1(3.5782337572542144, 2.411806254420634)
2(3.5328283559874127, 2.4714506345289475)
3(3.574915806878003, 2.087969844274367)
4(3.4637117664130024, 1.8035546342774778)
5(3.4688455610266042, 2.176686370574768)
6(3.5747851783476055, 2.5276117593106093)
7(3.4919292794303587, 2.326931254790753)
8(3.4371912907605564, 2.3372896519852655)
......
9(3.458757582712308, 1.859578874940672)
y_given_Xeq3p5 = ( Y | (abs(X - 3.5) < 0.1) ).sim(10000)

plt.figure()
y_given_Xeq3p5.plot()
plt.show()

The true conditional distribution of \(Y\) given \(X=3.5\) is Uniform(1.75, 2.5). We see that the approximate distribution is a little rough at the boundaries (near 1.75 and 2.5), because of the “close to 3.5” approximation. But we would not be able to simulate the conditional distribution given \(X=3.5\) without more sophisticated methods.

5.28 Solution to Exercise 3.28

sigma = 1

X, Y = RV(Normal(0, sigma) ** 2)

R = sqrt(X ** 2 + Y ** 2)

plt.figure()
(X & Y).sim(1000).plot()
plt.show()

plt.figure()
R.sim(10000).plot()
plt.show()

def dartboard_sim(sigma):
  X, Y = RV(Normal(0, sigma) ** 2)
  
  R = sqrt(X ** 2 + Y ** 2)
  
  r = R.sim(10000)
  
  return r.count_lt(1) / r.count(), r.count_gt(12) / r.count(), r.mean()
sigmas = np.arange(0.1, 5.1, 0.1)

results = np.array([dartboard_sim(sigma) for sigma in sigmas])

plt.figure()
plt.plot(sigmas, results[:, 0])
plt.xlabel('sigma')
plt.ylabel('Approximate P(X  < 1)')
plt.show()

plt.figure()
plt.plot(sigmas, results[:, 1])
plt.xlabel('sigma')
plt.ylabel('Approximate P(X  > 12)')
plt.show()

plt.figure()
plt.plot(sigmas, results[:, 2])
plt.xlabel('sigma')
plt.ylabel('Approximate E(X)')
plt.show()